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The most important characteristics of the fragmentation of heterogeneous solids is that the mass 
(size) distribution of pieces is described by a power law functional form. The exponent of the 
distribution displays a high degree of universality depending mainly on the dimensionality and on the 
brittle-ductile mechanical response of the system. Recently, experiments and computer simulations 
have reported an energy dependence of the exponent increasing with the imparted energy. These 
novel findings question the phase transition picture of fragmentation phenomena, and have also 
practical importance for industrial applications. Based on large scale computer simulations here we 
uncover a robust mechanism which leads to the emergence of energy dependence in fragmentation 
processes resolving controversial issues on the problem: studying the impact induced breakup of 
plate-like objects with varying thickness in three dimensions we show that energy dependence occurs 
when a lower dimensional fragmenting object is embedded into a higher dimensional space. The 
reason is an underlying transition between two distinct fragmentation mechanisms controlled by the 
impact velocity at low plate thicknesses, while it is hindered for three-dimensional bulk systems. The 
mass distributions of the subsets of fragments dominated by the two cracking mechanisms proved 
to have an astonishing robustness at all plate thicknesses, which implies that the non-universality of 
the complete mass distribution is the consequence of blending the contributions of universal partial 
processes. 
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I. INTRODUCTION 

Fragmentation into numerous pieces occurs when a 
large amount of energy is imparted to a solid within a 
short time [l| . Impact induced fragmentation of hetero¬ 
geneous materials is abundant in nature having also a 
high importance for industrial applications especially in 
mining and ore processing [2|-[8| . During the past decades 
research on fragmentation mainly focused on the statis¬ 
tics of fragment masses (sizes) obtained by the breakup 
of heterogeneous materials [il,S[i3. A large number of 
experimental [EIM and theoretical studies [TsI, [isl - l^ 
have confirmed that the mass distribution of fragments is 
described by a power law functional form. The exponent 
of the distribution was found to show a high degree of 
robustness, i.e. investigations revealed that the value of 
the exponent does not depend on the type of materials, 
amount of input energy and on the way the energy is 
imparted to the system until materials of a high degree 
of heterogeneity are fragmented [ESIii. The value of 
the exponent is mainl y d etermined b y th e dimensional¬ 
ity of the system [l3|, [III, [l8|, [l9|, [2l|, [HI, I25l4^ and by 
the brittle or ductile mechanical response of the material 
[HI- The universality of fragmenting has been shown to 
be the fingerprint of an underlying phase transition from 
the damaged to the fragmented phase of the breakup 

process [bl. 111. [13. [Tsj . 

Recently, experiments on the iuipact induced fragmen¬ 
tation of long th in g lass rods [HI, IHl and freely-hanging 
glass plates [sH, |32| revealed energy dependence of the 
mass distribution exponent, i.e. the exponent was found 
to increase with the imparted energy [HHHl • The impor¬ 


tance of these novel findings originates from the fact that, 
on the one hand, they question the universality and hence 
the phase transition interpretation of fragmentation phe¬ 
nomena, and on the other hand, they have consequences 
on the design of engineering technologies used for crush¬ 
ing in mining and ore processing M- Recent computer 
simulations have also provided an interesting counter ex¬ 
ample [HI- mass distributions of pieces obtained by the 
breakup of spherical bodies impacting against a hard wall 
have been found to get steeper with increasing impact 
velocity. However, it proved to be an apparent energy 
dependence which occurs solely due to the moving cutoff 
of the mass distributions, and hence, it can be trans¬ 
formed out by rescaling with the average fragment mass 
[HI- This study highlighted the importance of scaling 
and data collapse analysis when evaluating fragmenta¬ 
tion results of finite size systems. 

In order to resolve controversial issues on the energy 
dependence of the exponent of fragment mass distribu¬ 
tions, in the present paper we study the impact induced 
breakup of heterogeneous materials by large scale com¬ 
puter simulations. Our results demonstrate that energy 
dependence emerges when the fragmenting object is em¬ 
bedded in a higher dimensional space. Studying the frag¬ 
mentation of plate-like objects in three dimensions (3D) 
we show that energy dependence is obtained for low plate 
thicknesses, while it disappears for thick plates. The rea¬ 
son is that due to the interplay of the geometry of the 
sample and of the embedding space a transition takes 
place in the system between two fragmentation mech¬ 
anisms as the impact velocity is increased: At low ve¬ 
locities the crack structure is determined by the inter¬ 
ference of elastic waves resulting in an essentially two- 
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dimensional crack pattern with a regular structure. High 
velocity impact gradually excites cracking in the 3D bulk 
of the solid giving rise to a highly disordered crack struc¬ 
ture and a steeper decay of the mass distribution. In 
3D bulk samples the transition is hindered so that a 
unique exponent emerges. In spite of the observed non¬ 
universality of the complete mass distribution, identify¬ 
ing subsets of fragments dominated by different cracking 
mechanisms an astonishing universality of their mass dis¬ 
tributions is revealed at all plate thicknesses. Our study 
provides a robust scenario which leads to the energy de¬ 
pendence of mass distribution exponents in fragmenta¬ 
tion phenomena but it still underlines the importance of 
universality. 


II. DISCRETE ELEMENT MODEL OF 
FRAGMENTATION 

We investigate the fragmentation of plate-like objects 
induced by impact of a projectile in the framework of a 
discrete element model (DEM) developed recently [36|. 
The model has proven successful in reproducing key fea¬ 
tures of fracture processes of heterogeneous materials. 
Here we briefly summarize the main steps of the model 
construction based on Ref. [36|. Similar modelling ap¬ 
proaches have been also used in Refs. [^, [^. [ 3 ^. 

The sample is represented as a random packing of 
spherical particles which was generate d by sedimenting 
particles in a rectangular container [^, . The diame¬ 

ter d of the particles was sampled from a uniform distri¬ 
bution in a narrow range {d) — Ad/2 < d < {d) Ad/2, 

where (d) denotes the average diameter. The range Ad 
of diameter values was set as Ad/ (d) = 0.05. In the sim¬ 
ulations plate-like samples were constructed with a rect¬ 
angular basis of side length L and height H. Simulations 
were carried out with a fixed extension L = 30 varying 
the height of the sample H in the range H = 3 — 15 mea¬ 
sured in units of the average particle diameter (d). The 
total number of particles in the samples falls between 
5000 {H/ (d) = 3) and 25000 {H/ (d) = 15). 

In the model cohesive interaction of particles is pro¬ 
vided by beams which connect the particles along the 
edges of a Delaunay triangulation of the initial particle 
positions. In 3D the total deformation of a beam is calcu¬ 
lated as the superposition of elongation, torsion, as well 
as bending and shearing. Crack formation is captured 
such that the beams, modeling cohesive forces between 
grains, can be broken according to a physical breaking 
rule, which takes into account the stretching and bend¬ 
ing of contacts 

/^ inax(e^,0j) ^ ^ 
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Here Sij denotes the axial strain of the beam between 
particles i and j, while 0^, and Qj are the bending an¬ 
gles of the beam ends. The parameters Sth and Qth con¬ 
trol the relative importance of the two breaking modes 


FIG. 1. (Color online) Geometrical setup of the simulations: 
rectangular samples of a square shaped basis were considered 
in such a way that the side length L/ (d) = 30 of the square 
was fixed and the height H/ (d) of the sample was varied 
from 3 to 15. Here a sample is presented for id/ (d) = 11. 
A particle in the middle of the front side of the sample was 
selected, together with its neighbors it got an initial velocity 
pointing inward the sample. The cylinders connecting the 
particles represent beams and the white arrow indicates the 
direction of the impact velocity. The inset shows a closer view 
on a small segment of the sample. 


[ 21 I, [ 2 ^, [ 25 I, [sil, [13 • In the model there is only structural 
disorder present, i.e. the breaking thresholds are constant 
eth = 0.002 and &th = 2^, however, the physical prop¬ 
erties of beams such as length, cross section, and elastic 
moduli, are determined by the random particle packing. 
At the broken beams along the surface of the spheres 
cracks are generated inside the solid and as a result of 
the successive beam breaking the solid falls apart. The 
interaction of those particles which are not connected by 
beams, e.g. because the beam has been broken, is de¬ 
scribed by the Hertz contact law [38[ . The Hertz contact 
ensures that force can be transmitted through crack faces 
when they are pressed against each other. The fragments 
are defined as sets of discrete particles connected by the 
remaining intact beams. The time evolution of the frag¬ 
menting solid is obtained by solvit the equations of mo¬ 
tion of the individual particles [HI, [39| until the entire 
system relaxes meaning that no beam breaking occurs 
during one thousand consecutive time steps and there is 
no energy stored in deformation. For more details of the 
model construction and parameter settings see Refs. [H| . 

Impact loading was performed in such a way that a 
single surface particle was selected in the middle of one 
of the side walls of the sample. Together with its con¬ 
tacting neighbors it got an initial velocity vq pointing 
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FIG. 2. (Color online) Time evolution of a fragmenting plate 
of thickness H/ (d) =3 embedded in the three-dimensional 
space. The impact velocity v^jc — 0.2 is slightly above the 
critical point of fragmentation. Compressed beams are green 
while the stretched ones are red so that the propagation of a 
shock wave can be observed. Particles are colored according 
to the fragment they belong to while the color (grey) of frag¬ 
ments is randomly selected from a palette, (d) presents the 
hnal state of the system where the sample is reassembled by 
placing the particles back to their original position. 


towards the center of mass of the body. This is equiva¬ 
lent to an experimental setup where the impactor does 
not penetrate the target but it is sto ppe d after hitting 
the target surface as e.g. in Refs. [^, |^. The geomet¬ 
rical setup of the simulations and the loading condition 
of impact is illustrated in Fig. [H For the smallest plate 
thickness H/ (d) = 3 the impact site practically spans the 
cross section of the samples while for higher thicknesses 
the loading condition gets close to a point-like impact. 
Computer simulations were performed to determine the 
sound speed c of the model material. In the presentation 
of the results lengths and velocities are made dimension¬ 
less by dividing them with the average particle diameter 
(d) and with the sound speed c, respectively. 


III. DAMAGE-FRAGMENTATION 
TRANSITION 

In order to investigate how the overall geometry of the 
system affects the outcomes of the breakup process for 
each plate thickness H we carried out simulations vary¬ 
ing the impact velocity in a broad range. To accumulate 
statistics simulations were repeated 2000 times for each 
parameter set with different realizations of the structural 
disorder. As a representative example Fig. [2] presents the 
time evolution of a plate of thickness H/ {d) =3 gener¬ 
ated by an impact with initial velocity vq/c = 0.2. The 


sample breaks into a large number of pieces due to the 
shock wave generated by the impact resulting in a broad 
distribution of fragment sizes. The figure also shows the 
final reassembled body where fragments can be easily 
identified. Of course, the degree of breakup strongly de¬ 
pends on the value of the impact velocity vq: at low vq 
the sample just gets damaged around the impact site, 
i.e. some beams break and small fragments comprising a 
few particles are ejected but the main part of the body 
retains its integrity. To achieve complete fragmentation, 
where even the largest fragment is significantly smaller 
than the original body, the impact velocity has to exceed 
a critical value Vc. To quantify the degree of breakup 
we determined the average mass of fragments {M 2 1 Mi) 
as the ratio of the second M 2 and first Mi moments of 
fragment masses. The kth. moment M^ of the fragment 
mass in a single simulation is defined as 

= ( 2 ) 

i 

where rrii denotes the mass of single pieces, while rrimax 
is the largest fragment mass. The sum runs over all 
fragments. The ratio of the two moments M 2 /Mi was 
determined in single simulations and then it was aver¬ 
aged over fragmentation events at each impact velocity 
vq. The inset of Fig. [3] shows that gradually increas¬ 
ing vd the average fragment mass increases due to the 
creation of larger fragments. Since the largest fragment 
is always removed from the moments in Eq. (|2]), the de¬ 
creasing branch of (M2/Mi) is caused by the absence of a 
dominating piece. Hence, the position of the sharp max¬ 
imum can be identified with the critical value Vc of the 
impact velocity where complete breakup occurs. The re¬ 
sult demonstrates that the system undergoes a transition 
as the impact velocity is varied from the damage phase 
('T’o < '^c) characterized by the presence of a dominating 
piece, to the fragmentation phase (i;o > Vc) where no 
major fragment prevails. The existence of the damage- 
fragmentation transition has been verified for various 
types of systems both W experiments |3|-l8|, [l^ and com¬ 
puter simulations [l5|, |l8|, |^, [2314^ . 

In the inset of Fig. [3] the critical velocity Vc is an in¬ 
creasing function of H because the total mass of the sam¬ 
ple Mtot increases with the plate thickness Mtot oc L^H. 
The main panel of Fig. [3] presents that rescaling vq and 
(M2/Ml) of the inset with appropriate powers of the 
plate thickness H the results obtained at different thick¬ 
nesses can be collapsed on a master curve. The high 
quality collapse implies the validity of the scaling form 

(M2/M1) {vo, H) = ( 3 ) 

where (j){x) denotes the scaling function. The exponents 
were obtained numerically as a = 0.2(3) and p = 0.5(2) 
giving the best collapse in Fig. [3l It follows from Eq. m 
that the critical impact velocity Vc increases as a power 
law of the plate thickness 

VcOcH^. 


(4) 
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FIG. 3. (Color online) Inset: Average mass of fragments 
{M 2 1 Ml) as a function of the impact velocity vq for all plate 
thicknesses H considered. The main panel presents the same 
data rescaled with appropriate powers of the plate thickness 
H to obtain collapse of the curves. 


In the following we focus on the probability distribution 
of the mass of fragments p{m) in the fragmented phase 
to understand how it evolves with the impact velocity at 
different plate thicknesses. 



FIG. 4. Mass distribution of fragments at different impact 
velocities for four different plate thicknesses H/ (d): (a) 3, {b) 
5, (c) 11, (d) 15. In the damage phase the distributions are 
composed of two distinct parts, i.e. the large residues form a 
peak at m/Mtot ~ 1, while the small fragments have a rapidly 
decreasing distribution. The two regimes are separated by a 
gap which gradually disappears as the critical velocity Vc is 
approached from below. The dashed straight lines represent 
power laws of exponents (a) 1.7 and 2.4, (b) 1.7, (c) 1.9, and 
(d) 1.9. 


IV. MASS DISTRIBUTION OF FRAGMENTS 

Figure m presents the fragment mass distribution p{m) 
for several values of the plate thickness H at different 
impact velocities. A generic feature of the distributions 
is that in the damage phase {vq < Vc) p{m) is com¬ 
posed of two distinct parts: due to the presence of a big 
residue, a peak of the mass distribution is formed close to 
m/Mtot ~ 1 while the distribution of small pieces has a 
rapidly decreasing functional form. The two regimes are 
separated by a gap which gradually disappears as the 
critical impact velocity is approached from below. It can 
be observed that in the fragmented phase (i;o > '^c) small 
sized fragments have a power law mass distribution 

p{m)(xm~^, (5) 

which is followed by a cutoff regime. The power law over 
a broad range first occurs at the critical point. 

The important feature of our results presented in Fig. 
HJis that for low plate thicknesses H we observe a gradual 
increase of the mass distribution exponent from r = 1.7 
obtained at the critical point Veto t = 2A reached in the 
limit of high vq values (see Fig. SJ^a)). However, increas¬ 
ing the plate thickness in Figs. 01^6, c, d) this dependence 
of r on the impact velocity gradually disappears and for 
high plate thicknesses H/ {d) >11 only a single value of 


the exponent r = 1.9 remains. For increasing vq the mass 
of the largest fragment must decrease which may result 
in an apparent increase of the exponent simply due to the 
shifting cutoff of the distributions. However, contrary to 
Ref. our analysis showed that the change of r in Fig. 
[4|is the real behavior of p{m), it cannot be transformed 
out by rescaling with the average fragment mass. 

A. Fragmentation mechanisms 

Our computer simulations revealed that the observed 
dependence of r on the impact velocity is caused by 
a transition between two distinct fragmentation mech¬ 
anisms, which emerges due to the interplay of the geom¬ 
etry of the sample and of the dimensionality of the em¬ 
bedding space. It can be seen in Fig. [2] that immediately 
after impact the specimen gets damaged in the vicinity of 
the impact site, i.e. in a small volume starting from the 
surface all beams get broken and single particles (powder 
in the model) are ejected from the specimen. The impact 
loading generates a shock wave which gets gradually at¬ 
tenuated by the breaking of beams and by the expansion 
over a larger volume. At sufficiently high vq the com¬ 
pression wave can pass through the sample and reflects 
back with opposite phase as a tensile wave at the free 
boundaries of the rectangular specimen freely evolving 
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FIG. 5. (Color online) Mass distribution of fragments p{m) 
obtained at the critical point of fragmentation Vc for the small¬ 
est thickness H/ {d) = 3. Fragments giving dominating con¬ 
tribution to different ranges of p{m) originate from well de¬ 
fined spatial regions of the sample. These regions are high¬ 
lighted in a single sample in the inset using the same colors 
(characters) as for the corresponding ranges of the mass distri¬ 
bution. For comparison we also present the mass distribution 
at the highest impact velocity i^o/c = 0.5 where the power law 
regime is significantly steeper. The slopes of the two straight 
lines are r — 1.7 and r = 2.4. 


in the three dimensional embedding space. In this veloc¬ 
ity range the stress field which gives rise to the formation 
of extended cracks and final breakup is determined by the 
interference pattern of compression and reflected tensile 
waves. This mechanism has the consequence that for the 
limit of thin plates H L the breakup of the specimen at 
the critical impact velocity is caused by a relatively regu¬ 
lar crack pattern which is essentially two-dimensional. To 
demonstrate how this fragmentation mechanism works 
Fig. [5] presents the mass distribution p{m) of pieces at 
the critical impact velocity for H/ (d) = 3 together with 
a sample in the inset where particles of different spatial 
regions are highlighted by different colors. It can be ob¬ 
served that p{m) has a power law behavior over a broad 
range of fragment masses, however, it is decorated by dis¬ 
tinct maxima. Detailed analysis revealed that the emer¬ 
gence of the maxima is the fingerprint of the regularity 
of the two-dimensional crack structure of the plate-like 
object, i.e. fragments giving dominating contribution to 
a maximum always emerge in the same spatial region of 
the specimen. To make it clear in Fig. [5] we assigned 
colors (characters) to the maxima of p{m) such that in 
the inset the same colors (characters) are used for the 
particles of the spatial regions where the corresponding 
fragments originate from. It can be seen that the smallest 
fragments comprising a few particles {grey (a)) are gen¬ 


erated in the destructed zone close to the impact point, 
while larger fragments are formed deeper in the sample 
{cyan (c)). The interference of elastic waves generates a 
highly stretched zone along the surface of the plate which 
gives rise to the detachment of surface fragments both on 
the left and right sides {green (d)) and on the back side 
{blue (e)). The two corners of the front side of the sam¬ 
ple {red (b)) result in a slight local maximum of p{m) 
between the grey (a) and cyan (c) regions. The largest 
fragment {yellow (f)) controlling the cutoff of the distri¬ 
bution is created inside the specimen close to the back 
side with a shape elongated perpendicular to the direc¬ 
tion of impact. The two maxima of the yellow region of 
p{m) are caused by the fact that the yellow region breaks 
into two major pieces with a high probability. 

As the impact velocity gets high enough the overall 
two dimensional character of the crack pattern disap¬ 
pears and most of the cracks are created in the three- 
dimensional bulk of the material. This second fragmen¬ 
tation mechanism gradually becomes dominating with in¬ 
creasing vq. As a consequence, in Fig. m^a) the cutoff of 
p{m) shifts toward smaller m and the fraction of large 
fragments decreases in the mass distribution giving rise 
to a higher value of the power law exponent r. It follows 
from the above arguments that the dependence of the 
mass distribution exponent r on the impact velocity vq 
is caused by the gradual crossover from the planar two- 
dimensional to the three-dimensional bulk crack struc¬ 
ture. At the critical velocity the 2D character dominates 
while in the limit of high vq the crack pattern is com¬ 
pletely three-dimensional. In the intermediate velocity 
range both mechanisms are present so that the observed 
mass distribution is a blend of their contributions. The 
crossover is gradual in the sense that in the highly de¬ 
stroyed zone around the impact site the crack pattern is 
three-dimensional already at the lowest impact velocities 
which then spreads over the sample as vq increases. To 
have a clear view on the two limits of p{m) with different 
exponents of the power law regimes in Fig. [5] the mass 
distribution is also presented for the highest impact ve¬ 
locity voje = 0.5 we considered. The two straight lines 
of the figure represent power laws of exponents r = 1.7 
and T = 2.4. 


B. Superposition of subsets of fragments 

The changing exponent of the mass distribution is the 
consequence of the gradual crossover of the crack struc¬ 
ture in the sample as the impact velocity increases. The 
two- and three-dimensional crack structures favor frag¬ 
ment formation in different spatial regions of the sam¬ 
ple having also different extensions. Hence, in order to 
understand how the crossover emerges we identified sets 
of fragments according to their position in the sample 
and analyzed how their contributions to the complete 
mass distributions evolve with the impact velocity and 
plate thickness. The key feature of fragments is whether 




FIG. 6. (Color online) Identification of fragment subsets 
based on the bounding box of fragments and of the com¬ 
plete sample. Two plates are shown with different thicknesses 
and impact velocities (a) H/ (d) = 5 and v^jc — 0.23, (b) 
H/ {d) = 11 and vq/c = 0.3. For each subset the largest 
fragment is highlighted with different colors: light blue (light 
grey), green (medium grey), and red (dark grey) stand for 
the spanning, surface, and bulk fragments, respectively. The 
bounding boxes are indicated by the wire frames. 


they are created in the bulk or on the surface of the 
initial body. Since the sample surface is rather irregu¬ 
lar, to identify the position of fragments we construct 
their bounding box and compare the location of its cor¬ 
ner points to the bounding box of the original sample. 
For this purpose in the final state of the fragmentation 
process we reassembled the sample, as it is illustrated 
in Fig. [2j(i), and determined the bounding box of each 
fragment aligning the sides of the box with the edges of 
the sample. Based on the position and extension of the 
bounding box three types of fragments are distinguished: 

• Bulk fragment if the corners of the bounding box 
all lie inside the sample, i.e. their distance from the 
surface of the bounding box of the sample is greater 
than a threshold distance 0.2(d), the fragment is 
considered to be a bulk fragment. 

• Surface fragment if any of the corners, but not all 
of them, are within the threshold distance to the 
surface of the sample’s bounding box, the fragment 
is called surface fragment. 




FIG. 7. (Golor online) Mass distribution of subsets of surface, 
spanning, and bulk fragments together with the complete dis¬ 
tribution of all fragments for two plate thicknesses id/ (d) = 5 
and 15 in the upper and lower rows, respectively. In both 
cases the results are presented for two values of the impact 
velocity vq slightly above the corresponding Vc (left column) 
and for the limit of high speed impact (right column). 


• Spanning fragment those fragments which span 
the sample at least in one direction, are called span¬ 
ning fragments. For spanning fragments all the cor¬ 
ners of the bounding box lie in the vicinity of the 
sample surface. 

Figure [6] demonstrates the identification of the three sub¬ 
sets of fragments in two plates of different thicknesses. 
For each subset a single fragment is highlighted together 
with its bounding box. 

The spanning fragments are typically formed by cracks, 
which connect two opposite sides of the sample. Such 
cracks emerge due to the global interference pattern of 
elastic waves. In thin plates below and at the critical im¬ 
pact velocity Vc most of the mass is comprised in span¬ 
ning fragments; bulk and surface pieces can mainly be 
formed around the destroyed zone at the impact site. 
Above Vc the fraction of spanning fragments rapidly de¬ 
creases, however, they always have the largest mass so 
that the spanning fragments determine the cutoff of the 
mass distribution p{m) at any vq. Figure [7] presents the 
mass distribution of the three subsets of fragments for 
a plate of thickness iT/ (d) = 5 at the critical impact 
velocity together with the complete distribution. Note 
that the partial distributions are normalized such that 
their integral is equal to their fraction in the complete 
set of fragments. The Figure clearly demonstrates that 
the cutoff and the large mass regime of p{m) is controlled 
by fragments which span the sample in the direction per¬ 
pendicular to the plate. 











In thin plates surface fragments are generated from 
spanning fragments by cracks formed in the bulk of the 
specimen segmenting the spanning cracks. Due to this 
geometric constraint, in Fig. El^a) the low mass regime 
of p{m) is dominated by surface fragments, while bulk 
pieces have only a minor contribution. It follows that for 
thin plates the overall power law character of the com¬ 
plete mass distribution p{m) originates from a mainly 
two dimensional crack pattern. The power law regime of 
the complete mass distribution covers a broad range of 
fragment masses due to the spanning fragments and their 
daughter pieces on the surface of the sample. The expo¬ 
nent r of the power law regime was obtained by fitting 
r = 1.7 ± 0.05 (compare also to Fig. |4)). 

Fig. ^b) presents the same mass distributions at the 
limit of high impact velocities. Here the range of the 
power law regime of the complete mass distribution p{m) 
gets reduced and the exponent increases to a high value 
r = 2.4T0.07 (compare also to Fig.|4]). As the impact ve¬ 
locity increases bulk cracking gets activated which leads 
to a three-dimensional crack pattern with a highly disor¬ 
dered structure. Consequently, the difference of surface 
and bulk fragments disappears, they both have the same 
mass distribution with a high exponent in Fig. [3^6). Still 
the spanning fragments have the largest mass but they 
are only the remainings of the detached pieces. 

Fig. [Tfc, d) present the corresponding results for 3D 
bulk bodies with H/ {d) = 15. Compared to the plate¬ 
like objects of Fig. [Tl^a, b) one can observe that the span¬ 
ning pieces do not have a dominant role, at low velocities 
they are formed by detachment, while at high velocities 
they are just the corners of the sample. Spanning frag¬ 
ments give rise to distinct humps of the distribution at 
high m values, however, also surface fragments have con¬ 
tribution to the cutoff of p{m). Fig. |4] demonstrated that 
for bulk bodies the fragment mass distribution exponent 
practically does not depend on the impact velocity, it has 
a unique value r = 1.9 ± 0.05. According to Figs. [3^c, d) 
the reason of the constant exponent is that the relative 
fraction of surface and bulk pieces does not depend on 
the impact velocity. 

C. Universality of partial mass distributions 

A very interesting outcome of our study is that in spite 
of the velocity dependence of the exponent observed for 
the complete mass distribution, the partial distributions 
of the subsets of bulk and surface fragments exhibit a 
high degree of universality. Figure [ 8 ] presents the scal¬ 
ing plot of the mass distributions of surface and bulk 
fragments obtained at different impact velocities for two 
plate thicknesses. It can be observed that rescaling the 
distributions with appropriate powers of the impact ve¬ 
locity good quality data collapse is obtained in all cases. 
This scaling analysis demonstrates that the partial dis¬ 
tributions obey the scaling law 

p{m) =v^^{mv^), ( 6 ) 



FIG. 8. (Color online) Data collapse analysis of the mass 
distribution of subsets of surface and bulk fragments for two 
plate thicknesses H/ {d) = 5 and 15 in the upper and lower 
rows, respectively. Rescaling the distributions with appropri¬ 
ate powers of the impact velocity above the critical point good 
quality data collapse is achieved. The bold lines represent fits 
of the master curves with Eq. 0. 


where the exponent 7 depends on the plate thickness H. 
Note that due to the normalization of the distributions 
the same exponent 7 has to be used along the horizontal 
and vertical axis. In Fig. [ 8 ] best collapse was obtained 
with the exponents H/ (d) = 5: 7 = 1.3 (surface) and 
7 = 0.2 (bulk); H/ (d) = 15: 7 = 0.7 (surface) and 
7 = 0.25 (bulk). The scaling function T(x) was fitted 
with the functional form 

T(x) oc x~^ exp (— (x/x*)^), (7) 

where the exponent k and the characteristic scale x* only 
control the shape of the cutoff. The most remarkable fea¬ 
ture of the results is that best fits of the scaling function 
T is obtained with r = 1.7 and r = 2.4 for surface and 
bulk fragments respectively, for all thicknesses. This re¬ 
sult implies that the partial distributions exhibit univer¬ 
sality as it has been observed for a broad class of frag¬ 
mentation phenomena. The observed non-universality of 
the complete distribution of all fragments originates from 
the blending of the distributions of subsets of fragments 
whose contributions depend both on the impact velocity 
and on the plate thickness. 

V. DISCUSSION 

We investigated the impact induced breakup of het¬ 
erogeneous brittle materials in the framework of a three- 
dimensional discrete element model focusing on the mass 
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distribution of fragments. Based on large scale computer 
simulations we resolved recent debates on the universal¬ 
ity of the power law exponent of the mass distribution 
which is crucial both from scientific point of view and 
for industrial applications. Simulations were carried out 
to investigate the impact induced breakup of plate-like 
objects where both the thickness of the plate and the 
impact velocity were varied in a broad range. Our com¬ 
puter simulations revealed that for thin plates embedded 
in the three-dimensional space the power law exponent of 
the fragment mass distribution has a strong dependence 
on the impact velocity: power law is first obtained at 
the critical velocity of impact with an exponent r = 1.7 
which then gradually increases to r = 2.4 for high im¬ 
pact velocities. However, for 3D bulk samples a unique 
exponent is obtained r = 1.9, dependence on the impact 
velocity can only be pointed out for the cutoff of the dis¬ 
tributions. Note that the value of the exponent r = 1.7 
of p(m) falls close to the theoretical prediction of Ref. 
1^ . 1 ^ based on the branching-merging scenario of dy¬ 
namic cracks: if fragments are formed by the merging of 
branches of splitting unstable cracks a universal exponent 
of the fragment mass distribution r = {2D — 1)/D was 
predicted depending solely on the dimensionality D of the 
embedding space. For D = 3 the formula yields r = 5/3 
in the vicinity of our numerical result, although, in our 
case simulations did not reveal a branching-merging se¬ 
quence of cracks. The exponent r = 1.9 of 3D bulk sam¬ 
ples is consistent with other DEM results, e.g. the same 
exponent was obtained for the fragmentation of brittle 
spheres impacted against a hard wall [111, [2^ . 

The reason of the velocity dependent exponent is that 
due to the interplay of the geometry of the sample and 
of the dimensionality of the embedding space a crossover 
occurs between two different fragmentation mechanisms. 
In the vicinity of the critical impact velocity the crack 
pattern is essentially two-dimensional determined by the 
interference pattern of compressive and tensile waves gen¬ 
erated by the impact. This crack pattern has a high de¬ 
gree of regularity which gives rise to local maxima of the 
fragment mass distribution on the overall power law func¬ 
tional form. At increasing impact velocities bulk cracking 
gets activated so that the crack structure becomes three- 
dimensional with a high degree of randomness. 

A similar effect of the interference pattern of elastic 


waves has been observed for slender rods where fragmen¬ 
tation was induced by a hit at the free rod end. The mass 
distribution of pieces proved to have discrete humps at 
certain fractions of the buckling wave length simi¬ 
lar to what we obtained for plates. Studying the impact 
induced breakup of thin glass plates, in the experiments 
of Ref. [3l| an increase of the mass distribution expo¬ 
nent was reported with increasing impact velocity. The 
authors argued that the effect can be attributed to the in¬ 
crease of the fractal dimension of the crack pattern, i.e. 
as the crack structure gets more-and-more space filling 
the mass distribution exponent increases and approaches 
a limit value [3l|. Our results clarify the background 
of these experimental findings unveiling the underlying 
mechanism. 

Comparing the bounding box of fragments and that 
of the complete sample we decomposed the fragment en¬ 
semble into subsets of bulk, surface, and spanning pieces. 
The formation of these fragments is governed by different 
cracking mechanisms. Scaling analysis showed a striking 
universality of the mass distributions of bulk and sur¬ 
face fragments with strongly different exponents. The 
results imply that the velocity dependence of the expo¬ 
nent of the complete mass distribution at intermediate 
velocities vq is observed due to the mixing of the con¬ 
tributions of the subsets of fragments, where the mixing 
ratio depends on vq. Our results have the general conse¬ 
quence that energy dependence of the mass distribution 
exponent of fragmentation phenomena can be expected 
when a low dimensional object is embedded into a higher 
dimensional space allowing for the emergence of a tran¬ 
sition in the spatial structure of cracks generated by the 
initial shock wave. 
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